function[diplat] = get_dip_lat(z,f),
%double get_dip_lat(double z, double f){
%double incl,dipcolat,diplat;
incl = asin(z/f);
dipcolat = 0;
if (incl > 0), dipcolat = atan(2/tan(incl));end;
if (abs(incl) < 1.0e-10), dipcolat = 0;end;
  if (incl < 0), dipcolat = pi + atan(2/tan(incl));end;
  diplat = 180.0 * (pi/2.0 - dipcolat)/ pi;